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REPORT ON THE PERIOD MAY 1, 1992 — FEBRUARY 28, 1993 
A. Research Accomplishments 

I am pleased to report that we have greatly exceeded the original goals for the first year of this project. 

1. Mathematical Theory of Grain Alignment 

We have formulated a new mathematical technique for predicting grain alignment and polarization in molec- 
ular clouds. The new technique is completely general, in the sense that it can be used to calculate the 
efficiency of any alignment mechanism. It computes the “Rayleigh reduction factor,” i?, 1 by numerically 
integrating the “Langevin equation for Brownian rotation.” An algorithm for solving the Langevin equa- 
tion has been developed and thoroughly debugged via benchmark calculations on special cases where exact 
analytic solutions for R exist. The benchmark tests show that we are able to calculate R accurate to typi- 
cally 2% (3<r error) and that any systematic errors which may be present in our calculations correspond to 
spurious predictions of alignment which are well below the limits of detectability. We have implemented our 
algorithm on a parallel-processing machine with 2048 independent processors. Consequently, we now have 
the capacity to produce more results in one month than the combined results of all previously published 
work on grain alignment theory! This work has been submitted for publication. 

2. Super- Paramagnetic Alignment of Molecular Cloud Grains 

Before this project was funded, no calculations existed which could be used to compare near- and far-infrared 
polarimetric observations with the predictions of one theory — super-paramagnetic relaxation which is gen- 
erally considered to be a viable model for grain alignment in molecular clouds. We have carried out such 
calculations. The grains were modeled as refractory cores with ice mantles, where the core and mantle sur- 
faces were represented as confocal, oblate spheroids of arbitrary eccentricity. We included a proper treatment 
of the Barnett effect, rotational anelasticity, Larmor precession, gas-grain collisions, thermal evaporation of 
molecules from the mantle surface, super-paramagnetic absorption, and thermal fluctuations in the grain 
magnetization. We developed and incorporated a quantitative theory for the effects of thermal emission and 
infrared absorption on the alignment of nonspherical grains. We also carried out an exhaustive numerical 
study of the super-paramagnetic alignment scenario. Our calculations show that super-paramagnetic align- 
ment, in the absence of suprathermal rotation, is not consistent with the degree of alignment inferred by 
Lee and Draine (1985, ApJ, 290, 211) from an analysis of infrared extinction and polarimetry for the line of 
sight toward the BN object. This conclusion is independent of any assumptions about the uncertain grain 
magnetic properties. Our results were presented at the January, 1993 AAS meeting and a related paper 
will be submitted for publication before May 1. Calculations for another paper, on the analogous effects for 
prolate spheroids, are in progress. 

3. Theory of Grain Alignment by Ambipolar Diffusion 

The differential motion of ions and neutral particles in a partially-ionized plasma is called “ambipolar dif- 
fusion.” A charged grain in such a plasma will also drift through the neutrals, with a gas-grain drift speed, 
Vdt determined by the balance between the Lorentz and gas drag forces on the grain. Furthermore, if Vd is 
comparable to the gas thermal speed, then the grains will become aligned by a process known as “Gold's 
mechanism.” Gas-grain drift at thermal speeds is predicted in recent theoretical models of protostars (by 
Mouschovias and collaborators) and the neutral gas ring at the Galactic center (by Wardle and Konigl). It 
follows that grain alignment and polarization are implicit predictions of these models which can be used to 
test their validity. 

One of our first-year goals was the development of a quantitative theory for grain alignment by ambipolar 
diffusion. Such a theory has now been developed by us for spherical grains. An exact analytic solution 
has been found for the alignment of spheres subject to the effects of ambipolar diffusion plus paramagnetic 
or super-paramagnetic absorption and the other processes listed in item 1 above. This solution has been 


1 Thc Rayleigh reduction factor is a statistic of the grain angular momentum distribution which measures the efficiency of 
alignment. 
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used, in conjunction with a model for the grain optical properties, to generate theoretical predictions of the 
polarized, far-infrared emission from warm grains in a gas undergoing ambipolar diffusion. Our calculations 
show that detectable linear polarizations (« 1 - 5%) can be produced by this mechanism. A paper on 
the basic physics of these processes will be submitted before May 1. Our results will be applied in Year 2 
to calculate the far-infrared polarization which is predicted implicitly by the Wardle-Konigl models of the 
Galactic center gas ring. These predictions will be compared with far-infrared polarimetry of the ring by 
Hildebrand and collaborators, thereby testing the Wardle-Konigl models. 


4 . Publications 

Roberge, W.G., DeGraff, T.A., and Flaherty, J.E. 1993, “The Langevin Equation and its Application to 
Grain Alignment in Molecular Clouds,” submitted to ApJ. 

Roberge, W.G., DeGraff, T.A., and Flaherty, J.E. 1992, “Super-Paramagnetic Alignment of Core-Mantle 
Grains,” BAAS, vol. 24, no. 4, p. 1121 (abstract). 

DeGraff, T.A., Roberge, W.G., and Flaherty, J.E. 1993, “Magnetic Grain Alignment in Molecular Clouds,” 
to be submitted before May 1, 1993. 

Roberge, W.G., Hanany, S., and Messinger, D. 1993, “Grain Alignment by Ambipolar Diffusion,” to be 
submitted before May 1, 1993. 

B. Summary of Expenditures 


Item 

Budget 

Actual 

PI Summer Salary 

$5,954 

$0 

Grad Student Tuition 

4,417 

5,736 

Grad Student Academic Yr. Salary 

4,600 

9,200 

Grad Student Summer Salary 

3,220 

6,398 

Typing and Clerical Assistance 

350 

353 

Fringe Benefits 

1,324 

74 

Travel (domestic) 

1,200 

1,251 

Computer Usage 

2,000 

2,000 

Publication Costs/Page Charges 

1,500 

38 

Communications 

175 

163 

Indirect Costs 

10,260 

9,787 

Totals 

$35,000 

$35,000 


Differences between the actual and budgeted expenditures were due to the following: 

• The graduate student (Trudy DeGraff) had to be supported for a full year, not half a year as planned. 
To cover the increase in grad student support, the PI did not pay himself summer salary. This had no 
effect on the Pi’s level of commitment to the project and he did not receive summer salary from any 
other sources. 

• The PI supported 2 grad students during summer, 1992 instead of one. DeGraff worked on the 
mathematical theory of grain alignment and David Messinger worked on the theory of grain alignment 
by ambipolar diffusion. 

• The page charges for Paper 1 have not been billed at this date. 
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ABSTRACT 

We describe a new computational method for solving problems in grain alignment the- 
ory which uses numerical integration of the Langevin equation for Brownian rotation. The 
new method is completely general in the sense that it generates the solution to a Fokker- 
Planck equation with arbitrary diffusion coefficients. We derive accurate expressions for the 
diffusion coefficients of refractory grains with ice mantles, on the ad hoc hypothesis that the 
grains rotate with thermal kinetic energies. We include the effects of internal dissipation by 
Barnett relaxation or rotational anelasticity, Larmor precession, gas-grain collisions, thermal 
evaporation, and paramagnetic or super-paramagnetic absorption. We develop a quantita- 
tive theory for the effects of thermal emission and infrared absorption on the alignment of 
nonspherical grains and assess the relevance of these processes to polarization in molecular 
clouds. We document the accuracy of our computational method by comparing numerical 
solutions for the Rayleigh reduction factor of magnetically aligned grains with exact solutions 
for spheres and thin disks. 

Subject headings: polarization — interstellar: grains 
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1 INTRODUCTION 


Interstellar polarization is intimately related to the rotational dynamics of dust grains. 
Observational evidence shows unambiguously that, in both atomic and molecular regions of 
the interstellar medium, polarization by dust is caused by nonspherical grains whose axes 
are partially aligned with respect to B, the ambient magnetic field (e.g., Whittet 1992). 
Because the direction of each grain’s angular momentum vector, J , is tightly correlated with 
the direction of the grain’s principal axis of largest rotational inertia (Purcell 1979 [P79], see 
also §3.4), the polarimetric observations imply further that J must be partially aligned with 
respect to B. Thus, “grain alignment” means “angular momentum alignment.” It can be 
shown (Lee k Draine 1985) that the degree of linear or circular polarization from a cloud of 
dust depends on the grain angular momentum distribution only via the “Rayleigh reduction 

factor,” 

= | [(cos J ^>-i] (11) 

(Greenberg 1968), where is the angle between J and B and the angle brackets denote the 
average over an ensemble of identical grains. Thus, the Rayleigh reduction factor character 
izes the degree of alignment, with R varying from zero for an isotropic angular momentum 
distribution to unity for a distribution with J perfectly aligned parallel to B. All observa- 
tions of polarized extinction and emission by dust are presently consistent with R > 0, i.e., 
with the alignment of J preferentially parallel rather than normal to B (Hildebrand 1989; 
Whittet 1992). 

Grain alignment theory is a problem in statistical mechanics, the object of which is to 
predict the angular momentum distribution and thus R in terms of the various torques which 
act upon the grains. It is a remarkable fact that, more than forty years after the discovery 
of interstellar polarization, a convincing theory of the alignment mechanism remains to be 
established. While several physical processes are known which produce torques of the right 
type to align the grains (see the reviews by Hildebrand 1988a, b), it is unclear whether any 
of the processes which have been studied extensively can overcome the disalignment of J by 
random gas-grain collisions. For example, alignment by paramagnetic torques as described in 
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the classical Davis-Greenstein mechanism (Davis k Greenstein 1951 [DG]) can be ruled out 
unequivocally in H I clouds (Spitzer 1978; Whittet 1992) and at least some molecular clouds 
(Jones, Hyland, k Bailey 1984; Lee k Draine 1985). Similarly, it is difficult to reconcile the 
near-ubiquity of interstellar polarization with alignment either by gas-grain streaming (Gold 
1952; Purcell 1969), or by anisotropic radiation fields (Harwit 1970), both of which require 
rather special physical conditions. However, two variants of the classical DG mechanism may 
plausibly explain both the H I and molecular cloud observations (Hildebrand 1989): Magnetic 
alignment is greatly enhanced if the grains contain small ferromagnetic inclusions, thereby 
increasing the aligning torque to “super-paramagnetic” values (Spitzer k Tukey 1951; Jones 
k Spitzer 1967 (JS67); Purcell k Spitzer 1971 [PS71]; Duley 1978). Alternatively, or perhaps 
in conjunction with super-paramagnetism, the grains may rotate at suprathermal kinetic 
energies, thereby rendering J impervious to disalignment by gas-grain collisions (Purcell 
1979 [P79], Spitzer k McGlynn 1979; Johnson et al. 1981; Johnson 1982). At present, 
however, observational evidence regarding the super-paramagnetism and/or suprathermal 
rotation of interstellar grains is inconclusive (Hildebrand 1989). An unfortunate consequence 
of the ambiguous cause of grain alignment is that the relation between the polarization- to- 
extinction ratio and magnetic field strength in clouds is unknown. 

The wealth of polarimetric data which have been acquired from molecular cloud obser- 
vations during the past decade may provide important clues to the resolution of this problem. 
For example, polarization is observed in cold, quiescent globules where the dust-to-gas tem- 
perature ratio, Td/Tg, is known to be close to unity (e.g., Jones et al. 1984). It is possible 
that these observations are inconsistent with super-paramagnetic alignment, which vanishes 
in the limit Td/T 3 — * 1 unless the grains are also rotating suprathermally. Unfortunately, 
however, there is presently no way to assess this possibility, because there are no theoret- 
ical predictions on super-paramagnetic alignment in the absence of suprathermal rotation 
which can be compared with the observations. In particular, the last extensive calculations 
(PS71) preceded the discovery of physical processes (Dolginov k Mytrophanov 1976; P79) 
which affect the grain rotational dynamics in qualitative ways. Clues to the origin of grain 
alignment, are presumably also contained in other phenomena discovered via molecular cloud 
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polarimetry, such ss the relationship between the polarization-to-extinction ratio, P ( A, and 
the extinction in molecular clouds (Tamura et al. 1987) and the wavelength dependence of 
linear polarization both in the near-infrared continuum (Wilking et al. 1980, 1982; Martin 
& Whittet 1990) and in infrared grain spectral features (e.g., Aitken et al. 1985) to name 
just a few examples. 

This is the first paper in a two-part series on paramagnetic and super-paramagnetic 
grain al ignment in the absence of suprathermal rotation a scenario which we will refer to 
as “thermal alignment” for brevity. Our objective is to provide calculations which can be 
used to test this scenario, by comparing quantitative predictions with observations such as 
polarimetry of quiescent globules. In this paper we describe our numerical solution tech- 
niques (§2), discuss relevant physical processes in molecular clouds (§3), and describe some 
benchmark calculations which calibrate the accuracy of our calculations (§4). Our results 
are summarized in §5. In the second paper of this series (DeGraff, Roberge, & Flaherty 
1993, hereafter Paper II), we present an exhaustive exploration of the parameter space for 
thermal alignment and describe an unambiguous observational test. 

2 THE LANGEVIN EQUATION FOR BROWNIAN ROTATION 

In this section we introduce the Langevin equation as an alternative to the Fokker- 
Planck equation for calculations on grain alignment. For a more thorough description of the 
Langevin equation and its equivalence to the Fokker-Planck equation, the reader is referred 
to the excellent monographs by Chandrasekhar (1943) and Gardiner (1990). For discussions 
of numerical solution techniques, see Fox (1987), Gard (1988), and references therein. 

2.1 Equivalence of the Fokker-Planck and Langevin Equations 

The Rayleigh reduction factor is a statistic of the distribution function, for 

the grain rotational angular momenta at time t. In principle f can be found accurately by 
solving the Fokker-Planck (henceforth FP) equation, 

i=x ' y ' z (2 ' 1) 
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(JS67; Cugnon 1971), where sums are implied by repeated subscripts, 7, is the t-th cartesian 
component of J, and the mean torque, 

(A7i) = ^^j), i = x,y,z, (2.2) 

and “diffusion tensor,” 

(A7, AT,) = ( A ^ ~ 2 ) > = ( 2 - 3 ) 

are functions of J and t that can be found by evaluating the time averages denoted by angle 
brackets (see §3). The random variable A J is the incremental change in J which is caused 
by the total external torque during the time interval At. Following common practice, we 
refer to the mean torque and diffusion tensor collectively as the tt diffusion coefficients. 

Solving the FP equation directly is a formidable problem (JS67; Cugnon 1971, 1983, 
1985) because R is extremely sensitive to numerical errors in the calculated distribution 
function when the alignment is weak. For example, the linear polarization P ~ 1% which has 
been detected in the far-infrared and submillimeter emission from grains is consistent with 
#<0.1 in some regions (Dragovan 1986). Now suppose that we wish to calculate R — 0.1 
with modest accuracy, say 10%. Then it follows from eq. (1.1) that we must compute (cos j3) 
with an accuracy *** 1%. Solving the FP equation for by finite difference methods 

would obviously require a large number of grid points in angular momentum space to satisfy 
even this modest accuracy requirement. It seems likely that perturbation methods could 
be developed to handle the weak alignment limit, but then one would still need to perform 
a finite difference calculation to find R when the grains are strongly aligned, a case which 
is also of practical interest (Lee & Draine 1985). Monte Carlo simulations present a viable 
method for calculating R (Purcell 1969; PS71), but are much less efficient than the techniques 
described below. 

In this paper we describe a computational approach which works whether the alignment 
is strong or weak, yields R values accurate to typically ~ 1% with modest amounts of compu- 
tation, and is no more difficult to implement than a Runge-Kutta algorithm for integrating 
an ordinary differential equation. Our method relies upon the well-known mathematical 
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equivalence of the FP equation and the “Langevin equation,” 

dJi = Ai(J, t ) dt + t) dwj , i = x, y, z , (2.4) 

which is a stochastic differential equation for the random variable J(t). Here the angular 
impulse which is imparted to the grain during an infinitesimal time interval dt has been 
written as the sum of its mean value plus a fluctuating part corresponding, respectively, 
to the first and second terms on the right-hand side of eq. (2.4). As discussed below, a 
simple quadrature scheme can be used to integrate eq. (2.4) numerically as an initial value 
problem, yielding a pseudo- random, time-dependent variable, J(t), whose statistics converge 
to those of J(t) as discussed below. Moreover, it can be shown (e.g., Gardiner 1990) that 
the Langevin and FP equations describe the same random variable if one chooses A to be 
the mean torque, 

Ai = {AJi}, i=x,y,z , (2.5) 

and B to be the matrix square root of the diffusion tensor, 

(. BB t ) x] s {AJi AJ>), i,j = *,*,*, (2.6) 

where B T is the transpose of B. 

The statistical fluctuations in dJ are characterized by the random functions dwj in 
eq. (2.4). The requirement that J{t) should satisfy the FP equation uniquely determines 
{du?j} to be a vector of statistically independent “differential Wiener increments,” which 
are defined as follows. The Wiener process is a random variable, w(t), whose distribution 
function satisfies the one-dimensional FP equation 

df w _ 1 d 2 / w , 9 7 s 

dt ~ 2 dw 2 ' 

i.e., the diffusion equation. The finite Wiener increment is defined to be the change in w , 

Aui = w(t 0 + At) — w(t 0 ), (2.8) 

where to is some specified initial time, w(to) is assumed to be known with certainty, and 
At = t — to is an arbitrary but finite time interval. One may verify by direct substitution into 
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the one-dimensional FP equation that the Wiener increments have a Gaussian distribution, 



f Aw ( Aw) = (2* A t)~ l/2 exp [— ( Au>) 2 /2 At] , 

(2.9) 

with zero mean, 


(Aw) = 0, 

(2.10) 

and variance, 


> 

to 

II 

t> 

(2.11) 


a familiar result from the theory of random walks. The infinitesimal Wiener increment is 
simply the limit of Au> as A t — * dt. 


2.2 Numerical Integration Techniques 

Suppose that we wish to solve the Langevin equation for the components of J(t) on 
some time interval, [0, T], subject to the initial conditions 

Ji{ 0) = J.,o, i = x,y,z, (2-12) 

where may be known with certainty or may themselves be random variables. We 

partition the interval [0, T] into the uniform grid t n = nAt, n = 0, . . . , N, where N is an 
integer and At = T/N. In the simplest numerical method for integrating eq. (2.4), called 
the Euler-Maruyama scheme, we set Ji(t n ) 83 Ji,m where 

= i = x,y,z (2-13) 

and 

= Ji,n + Ai ( j(t„), t n ) At + Bij (J(t n ),tn) AtO,, i = X, V ,Z (2.14) 

(Gard 1988). The starting values are obtained by sampling from the appropriate distribu- 
tion if the initial conditions are random and the Wiener increments are obtained by sampling 

from the distribution given in eq. (2.9). Thus, the Euler-Maruyama scheme is a straightfor- 
ward generalization of Euler’s method for integrating nonrandom differential equations. The 
analogous Runge-Kutta schemes of higher order also exist, but they can be used only if the 
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diffusion tensor satisfies an additional symmetry condition and are not considered further 
here. 


Each pseudo-random variable J,(t) is a numerical simulation of the time-dependent 
“sample path” followed by as it traverses angular momentum space. If the diffusion 
coefficients satisfy certain smoothness and growth conditions, then Ji{t) converges to J x (t) 
in the mean square square sense, 

jirn^ ^ [./;(<) — 7, (<)] (2.15) 

(Gard 1988), where the angle brackets in eq. (2.15) denote expectation values at fixed t for a 
sequence of identical trial calculations, not time averages. The numerical solutions converge 
with the step size, At, at a rate such that the mean square error in J,, n is 

( [4. -•*((.))’) = 0( At”) (2' 16 ) 

where m = 2 for a single step (local error) and m = 1 for N steps (global error). The 
smoothness and growth conditions cited above are mild (Gard 1988) and are satisfied for all 
problems discussed in this and subsequent papers in this series (Paper II; Roberge & Hanany 
1993). 


For reasons discussed below (§4), we are interested only in the steady state angular 
momentum distribution. Thus we compute the Rayleigh reduction factor from the sequence 
{*/•,*», n = 0, . . . , in the obvious way, by evaluating cos 2 from the instantaneous value 
of the vectorial angular momentum at each time step and then taking the ensemble average, 
(cos 2 /?), to be the time average of cos 2 /? over a long sequence of consecutive steps. The time 
averaging begins at a time which is sufficiently large so that the grain has “forgotten” its 
initial conditions, which are therefore arbitrary. The total averaging time, T, is determined 
by the requirement that statistical fluctuations in the running average of (cos 2 /3) should be 
sufficiently small to guarantee the convergence of R to within a predefined tolerance. The 
exact value of T and other details of our calculations axe described in §4. 


Integrating the Langevin equation is similar to the Monte Carlo method in that both 
methods simulate the sample paths of J(t) and that the root-mean-square error in (cos 2 /3) 
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ia inversely proportional to the square root of the total integration or simulation time. 1 

1 This assertion is justified in §4 for the Langevin equation method. 

However, the Langevin equation has several significant advantages over the Monte Carlo 
technique. First, a computer code to integrate eq. (2.4) requires just a few lines of code to 
perform the time step in eq. (2.14) plus instructions to evaluate A, B, the Wiener incre- 
ments, and the running average of cos 2 (3. In contrast, the Monte Carlo method requires 
direct simulations of all relevant physical processes, each of which may require several steps. 
For example, simulating a gas-grain collision involves choosing the arrival time, thermal ve- 
locity, and collision site of the incident gas particle from the appropriate distributions and 
then simulating the stochastic dynamics of the collision. Second, the convergence theorems 
above guarantee that statistics such as R which are calculated numerically via the Langevin 
equation converge to the statistics of the exact distribution, /, and provide guidance for es- 
timating the dependence of statistical errors on A t and T . Finally, integrating the Langevin 
equation eliminates “coarse graining” procedures which are unavoidable in the Monte Carlo 
method. For example, the number of gas-grain collisions (typically ~ 10® during each corre- 
lation time of the random motion of a real grain) must be reduced to a manageable number 
by arbitrarily seeding up the masses of the gas particles in a Monte Carlo simulation. The 
only assumptions incorporated in the Langevin equation are those on which the FP equation 
itself is based— that J{t) is a Markovian random variable and has continuous sample paths. 

3 DIFFUSION COEFFICIENTS FOR THERMAL GRAIN ALIGNMENT 

In this section we calculate the diffusion coefficients for various physical processes on 
the ad hoc assumption that the grain rotational energies are of order kT g , where T g is the 
gas kinetic temperature. Analogous calculations for the alternative, suprathermal alignment 
scenario (P79) will be presented in a future paper. 
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3.1 Assumptions 


3.1.1 Grain Properties 

We adopt a two-component model in which each grain consists of a refractory core sur- 
rounded by a volatile ice mantle, both of which are at a common temperature, TV Because 
the grain composition is controversial, we give analytic expressions for the diffusion coeffi- 
cients wherein arbitrary grain models can be represented by appropriate choices of various 
material properties. Thus in the present paper we treat the i-th grain component (where 
* = c or m) as an unspecified, homogeneous solid with mass density complex magnetic 
susceptibility Xi = x( + * Xii complex dielectric function e^(A) at wavelength A. We 

assume that the imaginary component of x« is given by an expression of the form 

Xi= KiSl (3-1) 


(JS67), where Cl is the grain angular velocity and the quantity Ki sets the time scale for 
alignment by paramagnetic or super-paramagnetic torques (cf. §3.7). For a broad range of 
ordinary paramagnetic materials, 

2.5 x 10~ 12 Kelvin-s 0 , 

Ki « Y d 

(Purcell 1969; Spitzer 1978). For super- paramagnetic grains, Ki is larger by an uncertain 
factor < 10 3 — 10 s (JS67; Duley 1978) and the real part of x» is similarly enhanced over typical 
values for ordinary substances. We assume that the principle of superposition holds for the 


volume susceptibility, so that 


V C K C + V m K„ 


is the effective A"- value for the grain as a whole, where V c and V m are the core and mantle 
volumes, respectively, and V s V c + V m is the total grain volume. 


The shapes of interstellar grains are constrained by models of the wavelength-dependence 
of extinction and polarization in the vicinity of grain spectral features (Lee k Draine 1985; 
Draine 1988). The limited observational evidence which is presently available on molecular 
cloud grains (Aitken et ad. 1985; Lee k Draine 1985; Hildebrand 1988b) suggests that they 
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are “typically” oblate, with modest axis ratios ~ 2:1. In this paper we model the grains 
as oblate spheroids, there being no apparent reason for considering ellipsoids or other more 
complicated shapes. 2 Following the notation of Draine fe Lee (1984), we take the mantle 
2 That is, one can model observations of interstellar polarization successfully without resort- 
ing to more complicated shapes (Hildebrand 1988b). However the grain shape might have a 
dramatic effect on the dynamics of alignment. For example, Mathis (private communication) 
has noted that suprathermai rotation might be driven by gas-grain collisions and evapora- 
tion from the large-scale surface structures on a highly irregular grain. Consideration of 
this interesting possibility is deferred to a future paper on suprathermai spinup processes in 

molecular clouds. 

surface to be a spheroid with semiaxes a m parallel to the symmetry axis and b m > a m 
perpendicular to the symmetry axis. We assume for simplicity that the core-mantle interface 
is a spheroid confocal with the mantle surface and denote the core semi-axes by a c and b c . 
The eccentricity of the i-th component, 

ti = y/l - (a,-/fci) 2 » ( 3 - 4 ) 


may be different for i = c and t = m. The mass of our model grain is 




M d = -*t> m 

and its inertia tensor has components 

II ~ 

for rotation about the symmetry axis and 

Ixx ~ (l + a m/^m) 

for rotation about any axis perpendicular to the symmetry axis, where the factors 


©.* = 1 + 


— - 1 I (— ^ I — 1 
^ Pm / \bm / 


and 


©xx = 1 + 


p c f \ / b c \ 

\Pm J \bm / 


(i + <4/*U 


(3.5) 


(3.6) 


(3.7) 


(3.8) 


(3.9) 


vary from unity for homogeneous grains to pc! Pm f° r bare grains with no ice mantles. 
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3.1.2 Gas Properties 


Our model grains are assumed to be located in a cloud which has a constant, uniform 
number density, n, and magnetic field, B, and which is composed entirely of H 2 molecules 
with mass m. The gas particles are assumed to have a Maxwellian velocity distribution with 
kinetic temperature T g in a reference frame (the “gas frame” ) where the bulk velocity of 
the gas is zero. We solve the Langevin equation for the components of J in the gas frame 
under the assumption that the latter is an inertial reference frame; however, our calculations 
should also be accurate for accelerating plasmas whenever the time scale for acceleration is 
much longer than the time to establish the steady state angular momentum distribution. 
The latter is strictly less than the “gas damping time,” t ga , (cf. eq. [3.31]), which in turn is 
comparable to the time it takes a grain to collide with an amount of gas equal to its own 
mass. Numerical estimates of t gat and some other relevant time scales are given in Table 1. 
The center of mass (CM) velocity of the grain in the gas frame is assumed to satisfy v d < v th , 
where the gas thermal speed is defined to be 



Radiation pressure and other forces can cause a grain to drift through the gas at speeds 
Vd > in some environments, with significant consequences for the alignment process 
(Gold 1952; Purcell 1969; PS71). These effects will be studied in detail in a subsequent 
paper (Roberge & Hanany 1993) but are not considered here. 

3.2 Frames of Reference 

We solve the Langevin equation in the gas frame but in some cases it is easier to 
calculate the diffusion coefficients in a “body frame” whose axes are aligned instantaneously 
with the grain’s principal axes of inertia. Henceforth, vector and tensor components without 
superscripts will refer to the gas frame, whose cartesian basis vectors, (x, y, i), have z oriented 
along the magnetic field and x and y oriented arbitrarily in the plane normal to B. Vector 
and tensor components with a superscript “6” will refer to the body frame, whose basis 
vectors, (x 6 , y b , i 6 ) , have z b oriented parallel to the grain symmetry axis and the axes x b and 
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y b oriented arbitrarily in the equatorial plane. Henceforth we will refer to z b as the “major 
axis of inertia” because I b , > I b x . While the body frame axes are parallel to the principal 
axes of inertia, they are defined to be stationary in the gas frame and so do not corotate with 
the grain material. Thus the body frame is an inertial frame which is defined independently 
at each time by the instantaneous grain orientation. 

We specify the relative orientation of the gas and body frames in terms of the familiar 
Eulerian angles (Fig. 1). In setting one of the Eulerian angles equal to 0, we have anticipated 
the result (P79, §3.4) that J is closely aligned with the major axis of inertia at all times. 
The cartesian components of the mean torque and diffusion tensor transform according to 

(A Ji) = Gij (A J b ) (3-11) 

and 

(A Ji A Jj) = G im (A J b m A J*) G~). (3-12) 

where the transformation matrix is given in terms of the Eulerian angles by 


with 


and 


G = C D E 


C = 


cos <p — sin <p 0 
sin 4> cos <f> 0 

0 0 1 


\ 


D = 


V 


1 0 0 
0 cos 0 —sin 0 
0 sin/? cos 0 J 


^ cos tp — sin ip 0 X 
E = | sin ip cos ip 0 

0 0 1 


(3.13) 


(3.14) 


(3.15) 


(3.16) 


(Goldstein 1950). 


14 



3.3 Fluctuation-Dissipation Theorem 


There is a useful relationship between the diffusion coefficients for any physical process 
which tends, in the absence of other processes, to establish thermal equilibrium (e.g., Reichl 
1980). In particular, the mean torque about the major axis of inertia is given by 




(3.17) 


(cf. P79, eq. [1]), where the damping time, 


U 


amp 


2/;,tr 

((AJ,*) 1 )’ 


(3.18) 


depends on the diffusion tensor and the temperature, T , at which thermal equilibrium pre- 
vails. Analogous relations exist for rotation about the other principal axes, but are not 
relevant here for reasons discussed below. Of course the alignment of interstellar grains is 
inherently a nonequilibrium situation, due to the different temperatures which are generally 
associated with different processes. Nevertheless, one can use the Fluctuation-Dissipation 
Theorem to find the mean torque for any individual process, such as thermal emission (see 


§3.8.3), which is characterized by a unique temperature. 


3.4 Internal Dissipation of Rotational Energy 

This term refers collectively to processes first described by Purcell (P79) which align 
the major axis of inertia with J. According to Purcell’s analysis, the most efficient of these 
processes is associated with the Barnett effect, in which the transfer of angular momentum 
from bulk rotation to aligned spins or orbits produces a magnetization 

Af = x ft/l (3.1.9) 

in a grain rotating with angular velocity ft. Here 7 is the magnetogyric ratio of the aligned 
magnetic moments. Purcell observed that the magnetization in eq. (3.19) is identical to that 
which would be produced by applying a “Barnett equivalent” magnetizing field H Be = ft / 7 
and, further, that if J is not aligned along a principal axis of the grain then ft, and hence the 
effective applied field, are nonsteady in a reference frame fixed to the grain. Magnetization 
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of the grain material by such a nonsteady “field” leads to paramagnetic absorption and the 
dissipation of rotational energy as in the Davis-Greenstein mechanism, except that the roles 
which are played in the DG mechanism by B and J are played in Barnett relaxation by fi 
and the major axis of inertia, respectively. Further, because Barnett relaxation involves no 
external t orques, J is a conserved quantity. 3 The transformation of rotational kinetic energy 

3 This statement is not strictly true. The quantity which is conserved is obviously the 
total grain internal angular momentum, F, which includes not only rotation but also the 
angular momentum 5 associated with the magnetization. However, for a homogeneous, 
oblate spheroid with equatorial radius b m composed of a paramagnetic or super-paramagnetic 
substance with partially aligned electron spins, one can show that 

S/J = 5Xm^C 2 /2e 2 pm^ « 10- 4 Xm/>m!o*C-5> 


where e and m e are respectively the electron charge and mass and the notation p m , o means 
p m j 10° g cm -3 and similarly for the values of other quantities in Gaussian cgs units. Thus, 
even for super-paramagnetic grains, the ratio of spin-to-rotational angular momentum is only 
S/J< 1%. We will therefore ignore the distinction between J and F, an assumption which 
is also implicit in Pu rcell’s (P79) discussion of Barnett relaxation. 

into heat by Barnett relaxation drives the grain toward the state of smallest rotational energy 
consistent with a given value of J, i.e., toward rotation about the major axis of inertia (P79). 


For interstellar grains, the time scale, for Barnett relaxation is always many orders 
of magnitude longer than the grain rotational period (Table 1). For spheroidal grains the 
motion of fi in grain coordinates is therefore essentially that of a free symmetric top over 
times small compared to tbar- From this information one can easily compute the energy 
dissipation rate and from it the rate of change of the angle 9 between fi and the major axis 
of inertia. Applying Purcell’s analysis to the case of an oblate, spheroidal, core-mantle grain, 
we find that 


d9 

dt 


VKh\h - 1)J 2 


sin 9 cos 9, 


7 2 (4) 3 

where K is defined in eq. (3.3), V is the total grain volume, and 


(3.20) 


h = I b JlL 


(3.21) 
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We will define the Barnett alignment time to be 


tba 



(3.22) 


A numerical estimate for molecular cloud grains is given in Table 1. As Purcell pointed 
out, tfar is typically so short in comparison to the time scales for J to be changed by 
collisional or magnetic torques, that one may consider the axis of largest rotational inertia 
to be aligned with J at all times. This conclusion is reinforced if one includes the other 
dissipative mechanism discussed by Purcell, rotational anelasticity, although the latter is 
typically one to two orders of magnitude less efficient than Barnett relaxation (P79). We 
note that this is the only consequence of internal dissipation on the alignment processes 
considered here because, as discussed in the footnote above, Barnett relaxation and rotational 
anelasticity conserve J to an excellent degree. 


3.5 Larmor Precession 


Barnett magnetization endows a spheroid rotating about its major axis of inertia with 
a steady magnetic moment 

.,'T / 

(3.23) 


_ X'V 

t * iar ~ rb **’ 


where = {V c x' c + V'mXjJ/V'. The magnetic torque /i 6ar xB causes J to precess about B 
at the Larmor frequency, fh a r = 2 tt / f( ar , where the Larmor period is 


2*7 /« 
X'VB- 


(3.24) 


There is also a relatively small contribution to the magnetic moment from rotating surface 
charges. For an oblate spheroid with charge Q uniformly distributed over its surface, the 
neglected magnetic moment would be smaller than the Barnett moment by a factor 


? (3.25) 

M6or 8?r ex' {a m /b m )b m ’ 

where X is a n um erical factor which varies from 2/3 for spheres to 1/2 for thin disks. The 
ratio (3.25) is typically only 0(1O -5 ) and we will neglect n c h g . Representative values of the 
Larmor period are given in Table 1. 
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The Larmor period is always several orders of magnitude shorter than the time scales for 
J to be changed by the collisional, paramagnetic, and radiative torques discussed below. It 
is therefore appropriate to average the diffusion coefficients for these torques over a uniform 
distribution of the Larmor precession angle, <f>, in Fig. 1. This is the only effect of Larmor 
precession on our calculations: once the diffusion coefficients for the other torques have been 
Larmor-averaged, it becomes unnecessary to include the torque n baT X B explicitly in the 
Langevin equation. (Note that Larmor precession has no effect on internal dissipation, even 
though the Larmor period can be shorter than the Barnett alignment time scale. This is 
because rotation about the axis of largest rotational inertia is the state of minimum kinetic 
energy for a given value of J and the Larmor torque is conservative.) 

For each of the processes considered below, the diffusion tensor is diagonal in the body 
frame with cartesian components 

' <(. Mi?) 0 0 ^ 

(AJjAJj) = 0 ((AJi) 2 ) o , ( 3 - 26 ) 

, 0 0 ((AJ{) 2 ) j 

where the subscripts “|| B and “1,” respectively, denote body frame axes which are parallel 
and normal to- the symmetry axis. Transforming these components to the gas frame using 
eq. (3.12) and averaging over the phase angles <j> for Larmor precession and ip for rotation, 
we find the diffusion tensor in the gas frame. The latter is also diagonal with components 

((A J.)*} . i { (1 + cos’ f» ((AJi)’) + sin’ ft <(A •/{)’) } , (3.27) 

{(AJ,)’) = {(AJ I )’), (328) 

and 

((AJ,)’) = sinV((AJl)’) + cos’/j((A (3.29) 

3.6 Gas-Grain Collisions and Thermal Evaporation 

The diffusion coefficients for gas-grain collisions and thermal evaporation have been 
calculated in Appendixes A-C under the following assumptions: (1) We assume that every 
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colliding H 3 molecule sticks to the surface long enough (i.e., a few times the surface oscillation 
period for physisorbed Hj) so that its kinetic energy becomes thermalized at the grain 
temperature; (2) we assume that every sticking collision is followed by thermal evaporation 
and that the velocity distribution of evaporating molecules is determined by the principle 
of detailed balancing (Burke k Hollenbach 1983); (3) we assume that evaporation occurs 
uniformly over the grain surface, i.e., that there are no preferred evaporation sites; and (4) 
we set the gas-grain drift velocity to v d = 0 as discussed above. Assumptions (l)-(2) should 
be excellent approximations for H 2 colliding at thermal energies either with cold grains (with 
3 < Tj < 15 K) covered with H 2 or with warmer grains whose ice mantles are devoid of H 2 
(see the discussion in Tielens k Allamandola 1987). Assumption (3) is consistent with the 
thermal alignment hypothesis under consideration in this paper, but might be inappropriate 
for real grains if the latter are highly irregular (see the footnote in §3.1.1). Nevertheless, item 
(3) is plausible because any preferred evaporation sites, if present, will generally have lifetimes 
that are smaller than the time to establish suprathermal rotation (which is comparable to 
t gat ). For example, thermal diffusion of the H 2 molecules on a cold grain will rearrange the 
surface topography on a time scale, t 1 s, which is always negligible compared to t ga3 . For 
warmer grains, the corresponding time scale for H 2 0 diffusion on an H 2 0 surface is still less 
than t gat if T d > 30 K. See Tielens k Allamandola (1987) and references therein. 

Within the above approximations, the the mean torque for collisions plus evaporation 
has gas frame components 

(A4U. = -•*/<».. ( 33 °) 


where 

t = — ^ (3.31) 

" 4v^F nmt^r|,(e m ) V 

is the gas damping time. The diffusion tensor is diagonal in the body-frame with components 


and 


((AJ|f) 2 ) c+ev = ^-nrn 2 b 4 m v* h ^1 -I- r ll( e ”0 ( 3 - 32 ) 

<(*4) 3 > c+ev = (* + t) < 3 - 33 ) 
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The only mathematical approximation used in deriving expressions (3.30) (3.33) is the ne- 
glect of terms of order y/m/M d and higher (Appendix C). The grain shape enters only via 
the “eccentricity factors,” 

r„(e m ) = ^ {3 + 4(1- e 2 Jg(e m ) - e“ 2 [l - (1 - e 2 J 2 g(e m j \ } (3.34) 

and 

TxM = 4 { 7 - + (■ - 4.)jM + (1 - 24 .) [1 + [1 - U - 4fsW]|} . 

(3.35) 


where 



(3.36) 


and e m is the eccentricity of the mantle surface. These functions have the limiting values 


lim r.|(e m ) = lim T x (e m ) = 1 


(3.37) 


for spheres and 

lim rii(e m ) = lim Tx(e m ) = - ( 3 - 38 ) 

for infinitely thin disks. The eccentricity factors are plotted in Fig. 2. One can easily 
verify that the diffusion coefficients in eq. (3.30)-(3.33) satisfy the Fluctuation- Dissipation 
Theorem and reduce to the expressions derived by Jones & Spitzer (JS67) for spheres in the 
appropriate limit. 

We obtain the gas frame components of the diffusion tensor by substituting eq. (3.32)- 
(3.33) into eq. (3.27)-(3.29) with the result that 

((AJ,)’} c+jv = (l + y ) [(l+cos^Ix + rin’/H’,], 

and 

({AJ,) 2 ) c+ev = ^nm 2 b A m v* h (l + ^ [3in 2 ^r x + cos 2 ^^,]. 


(3.39) 

(3.40) 

( 3 -41) 
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3.7 Paramagnetic or Super-Paramagnetic Absorption 


Jones & Spitzer (JS67) derived the diffusion coefficients for paramagnetic or super- 
paramagnetic absorption in a spheroidal grain. The mean torque has gas frame components 


— ~Tijt ma gy * — x iVi 

(3.42) 

and 


(AJ,) mag = 0. 

(3.43) 

The magnetic damping time for a core-mantle grain is 


I h 

'■mag — gy q 2 ' 

(3.44) 

The diffusion tensor is diagonal with gas frame components 


((A Ji) 2 ) =2kT d VB 2 K, i = x,y , 

(3.45) 

and 


((A J,)*) = 0. 

\ / mag 

(3.46) 


These diffusion coefficients satisfy the Fluctuation-Dissipation Theorem by construction. 


3.8 Electric Dipole Emission and Absorption 

The an gular momentum of a grain is changed by the emission, absorption, and scat- 
tering of photons. Purcell t Spitzer (PS71) pointed out that radiative processes generally 
contribute to both the mean torque and diffusion tensor and estimated the diffusion co- 
efficients for spherical grains in equilibrium with a thermal radiation bath at temperature 
Trad — Td' In view of the growing body of observations on polarized thermal emission from 
dust, especially the warm dust with Td ~ 100 K near the Galactic center (Werner et al. 1988; 
Hildebrand et al. 1990, 1992), it is appropriate to give a quantitative description of these 
effects for nonspherical grains. We therefore calculate the diffusion coefficients for emission 
from and absorption by oblate, spheroidal, core-mantle grains and consider explicitly the 
nonequilibrium case where T ra d ^ TV The only assumption used in the following discussion 
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is that the grain-photon interaction is in the Rayleigh limit. This is a good approximation 
for bare silicate grains with radii b c < 1 /jm and graphite grains with radii b c < 0.1 /im 
if the photon wavelengths are A >10 /im (Draine & Lee 1984). The following discussion 
therefore applies mainly to grains that are shielded from direct sources of optical and UV 
photons. Grain alignment and disalignment by short-wavelength radiation (Harwit 1970), 
which might be important in photodissociation regions, is not considered here. 

3.8.1 Electric Dipole Cross Section of a Rotating Dielectric Spheroid 

In the electric dipole approximation, scattering is negligible and the gram absorption 
cross section at frequency u is given by a simple analytic expression (van de Hulst 1981). 
The cross section for a non-rotating grain interacting with a linearly polarized plane wave 
whose electric field vector, E , lies parallel to the j-th principal grain axis is 

cm-*? imK) i hi.-l. < X47 > 

where u> is the circular frequency of the wave and is the electric polarizability of the grain 
for E polarized along the j-th body frame axis. The polarizability tensor of a spheroidal 
core-mantle grain can be expressed in terms of the dielectric functions of the core and mantle 
materials, the grain eccentricity, find the core-to-mantle volume ratio. For details, see Draine 
& Lee (1984). 

There is a small correction to expression (3.47) when the grain rotates, which arises 
because the dielectric polarization of the grain depends on the frequency of E relative to the 
grain material. Consider an oblate spheroid rotating about its symmetry axis with angular 
velocity fl = J/I b zx . Let E Q be the electric field amplitude of an incident plane wave and let 
£ be the angle between the electromagnetic wavevector, k , and the grain rotation axis, z b . 
It is useful to take the two independent polarization states of the incident wave to be the 
states of positive and negative helicity, corresponding to left and right circular polarization, 
respectively. Define a set of rotating basis vectors, (x r ,y r ,i r ), which are fixed to the grain 
material with z T along the rotation axis and with x r and y r oriented arbitrarily in the grain s 
equatorial plane. One can easily show by a suitable transformation of coordinates that the 
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electric field components in the rotating frame attached to the grain are 

££(±) = ^E 0 { [cos C T 1] cos (u> + fi) t + [cosC± l]cos(u> - fi)t } , (3.48) 

El(±) = -\e 0 { [cos C T 1] sin (w + fi) t - [cos C ± 1] sin (u - fi) t } , (3.49) 

and 

E T Z {±) = -E 0 sin C cos u>t, (3.50) 

where the upper (lower) sign choice corresponds to the case where the incident wave has 
positive (negative) helicity. Thus in general the grain material “sees” three waves: a wave at 
the incident frequency, w, which is linearly polarized along the grain rotation axis and two 
waves at frequencies wifi which are circularly polarized in the plane normal to the rotation 
axis. If we assume the grain material to be a linear dielectric, then the cross section for the 
g um of these waves is just the sum of the cross sections for the three frequency components 
weighted by their respective squared amplitudes. It folows that 

C(C, ±) = i [cosC T l] 2 Cx{u + fi) + ^ [cosC ± l] 3 c±(w - fi) 

+ i sin 2 C Cj|(u), (3-51) 

is the electric dipole absorption cross section for a circularly polarized wave propagating at 
angle C with respect to the symmetry axis, where the upper and lower sign choices correspond 
to the positive and negative helicity states of the incident wave. 

3.8.2 Absorption in an Isotropic Radiation Field 

Expression (3.51) implies that a rotating grain absorbing photons from an isotropic 
radiation field will be slowed by rotational friction. The effect, which was first described by 
Purcell & Spitzer (PS71), is illustrated by the special case of a wave propagating parallel 
to the rotation axis (C = 0). According to eq. (3.51), the absorption cross section for such 
a wave is C(0, +) = C±(u> — fi) if the wave has positive helicity and C( 0, — ) = C±(u + fi) 
if it has negative helicity. Of course this is just what one expects in this special example: 
the incident wave has E circularly polarized in the equatorial plane of the grain with E 


23 


rotating in the same (opposite) sense as the grain if the helicity is positive (negative). If 
the el ectric dipole cross section increases with frequency* then the grain is more likely to 

4 At frequencies where this assumption is violated, the torque due to electric dipole absorp- 
tion tends to spin up a rotating grain. Suprathermal spinup will occur if (a) C x decreases 
with u> over a sufficiently large fraction of the frequency range where photoabsorption is 
important and (b) the photoabsorption damping time is small compared to the damping 
time scales for other processes. However, neither (a) nor (b) is satisfied by the model silicate 
grains considered below. 

absorb the negative helicity wave. Recalling that in the electric dipole limit the photon 
orbital angular momentum is negligible compared to its spin angular momentum and that 
the spin of a positive (negative) helicity photon is parallel (antiparallel) to k, we see that 
the cross section is larger for the absorption of photons which tend to decrease J. This is 
the rotational friction noted above. 

Now consider the mean torque, (AJ) abl , due to photoabsorption in an unpolarized, 
isotropic field. (The generalization to anisotropic radiation is obvious and will not be con- 
sidered explicitly here.) The contribution to (AJ) o6< from a pencil beam of photons with 
wavevectors in an infinitesimal solid angle dW centered about the unit vector k and frequen- 
cies in the infinitesimal interval du> centered about u has body frame components 

<i(AJ‘) ali = [C«,+)-C((,-)l kldWJu (3.52) 

where I w is the monochromatic specific intensity of the radiation in photons cm 2 s 1 sr s. 
The term in eq. (3.52) which depends on the grain cross sections can be expressed in terms 
of the electric dipole cross sections for linearly polarized light using eq. (3.51). The algebra 
is simplified if we note first that u in all cases of practical interest. For example, if we 
consider a homogeneous grain with rotational energy equal to the equipartition energy at 
the gas temperature and take a typical photon energy to be kT ra d , then 

- ~ 10-* T,f p-„f b-JL\ (l - el) - ,/4 T-J. (3.53) 

CJ x 

The ratio in eq. (3.53) wll therefore be much less than unity even if the grain rotates at a 
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highly suprathermal energy. Thus, to within an excellent approximation we may set 

Ci(u> ± Q) ss Cx(w) ± fl (3.54) 

in expression (3.51), where 

Ci = ( 3 . 55 ) 

ajjj 

Making use of approximation (3.54), we find that 

C«, + )-C(C,-) « — cos C- (3-56) 

Substituting expression (3.56) into eq. (3.52), writing out the body frame components of k 
and performing the integrals over solid angle, we find that 

= —J/tabil (3.57) 

where the photoabsorption damping time is 

U” (3 - 58) 

In setting the upper limit of the frequency integral above to infinity, we are assuming the 
radiation temperature to be sufficiently small so that the integrand is only significant at 
frequencies where the electric dipole approximation is valid. For a thermal radiation field 
with temperature T rad < 100 K, the fraction of photons with wavelengths A < 10 /xm is less 
than 10 -4 . 

The incrementad contribution to the diffusion tensor from the monochromatic pencil 
beam discussed above haus cartesian components 

d (A J?A J ) = J J„ [ C(C, +) + C(C, -) ] it dW <L> (3.59) 

in the body frame. Substituting expressions (3.51) for the cross sections into eq. (3.59), 
writing out the body frame components of k, and performing the integral over solid angle, 
we find that the diffusion tensor is diagonad with nonzero components 

( (A ^ /” [C||M + 4CxM] L du (3.60) 
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and 


(3.61) 


((A/i)’) =^f 1°° [2C||H + 3Ci(»)] /„i*. 

One can easily verify that in the appropriate limit expressions (3.60) and (3.61) reduce to 
the corresponding expression derived for spherical grains by Purcell & Spitzer (PS71, eq. 
60). 

3.8.3 Thermal Emission 


The diffusion tensor for thermal emission can be deduced by noting that, in the hypo- 
thetical situation where the grain is immersed in a blackbody field with T ra d = Td, the prin- 
ciple of detailed balancing requires that emission and absorption should contribute equally 
to the diffusion tensor. The components are therefore 

((AJ,;) ! ) jT [CiiM + dCaH] a,(H)*>. (3.62) 

and 

((Adi) 2 ) JT < 3 - 63 ) 

where 

B„(T d ) = ( ex P (Ku/kTd) - l]' 1 (3-64) 

is the specific intensity of a thermal radiation field with temperature Td. 

The mean torque contributed by thermal emission can be computed using the Fluctuation- 
Dissipation Theorem. Thus in thermal equilibrium at temperature T d , the damping time for 
rotation about the symmetry axis due to absorption plus emission is 



where we have used the fact that thermal emission contributes half of the total diffusion 
tensor in thermal equilibrium. But because the mean torques due to absorption and emission 
axe additive, the mean torque due to emission alone must be 




(3.66) 
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(3.67) 


where the damping time for thermal emission is related to the other time scales by 

i~l _ 1 _ 

r em “ l TE l abf 

We have plotted the photoabsorption and photoemission damping times in Fig. 3 and 
Fig. 4, respectively, for the special but illustrative case of bare silicate grains with no ice 
mantles. We used the “astronomical silicate” dielectric function of Draine (1987) to compute 
the electric dipole absorption cross sections and chose p c — 2.5 g cm 3 . We see from the 
figures that in clouds where t gat < lO 10 s, the time scales t ab , and are at least an order 
of magnitude larger than t ga , unless the grains have radii b c < 0.1 /zm and the radiation 
and dust temperatures are larger than 100 K. We conclude that photoabsorption and pho- 
toemission typically have a negligible influence on the alignment of large dielectric grains 
in warm clouds, but might be important in studies of polarized emission from very small, 
transiently-heated grains. 


4 THERMAL ALIGNMENT IN MOLECULAR CLOUDS 

To assess the accuracy of our numerical methods, we consider a special case— the 
magnetic alignment of spheres and thin disks— where our numerical results can be compared 
with exact analytic solutions for R. Below we include all of the processes discussed in §3 
except for photoabsorption and thermal emission. 


4.1 Magnetic Alignment of Core-Mantle Spheroids 


The number of independent physical parameters in the Langevin equation is minimized 
by choosing appropriate dimensionless variables. Thus we define the dimensionless angular 
momentum components in the gas frame by 


Ji = 


Ji 




i — x,y,z, 


(4.1) 


and measure time in units of the gas damping time. When the diffusion coefficients from §3 
are written in terms of these variables and the resulting expressions for collisions, evaporation, 
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and paramagnetic or super-paramagnetic absorption are added together, one finds that the 
total mean torque has gas frame components 


(A;,) = -(1 + S)j„ 


(4.2) 


and 


(&jz) — —jzi 

(4.3) 

where 

« = ™ VB1 (4.4) 

iy/Srnmvthb^Tw 

is the ratio t ga Jt mag of the gas and magnetic damping times. The nonzero gas frame com- 

ponents of the dimensionless diffusion tensor are 


((Aj.) > = 2 (l + jr) [(l + cos f>) r|i(em) + »« <»] + 2 (j 

r) 3. (4-5) 

((Ajv) 2 ) = ) > 

(4.6) 

and 

«^>’>=H) H'SK + c H- 

(4.7) 


The coefficients which appear in the Langevin equation, eq. 2.4, can be obtained from the 
diffusion coefficients above via eq. (2.5)-(2.6). 


According to eq. (4.2)-(4.7), the diffusion coefficients, and hence the Rayleigh reduction 
factor, depend only on three dimensionless parameters. These are 6 , Td/T g , and the eccen- 
tricity of the mantle surface (which determines the ratio Tx/ru). Our results are identical 
in this respect to the earlier work of Purcell and Spitzer (PS71), even though the present 
work differs from PS71 in that we include Barnett relaxation and Larmor precession. Notice, 
however, that 6 as defined in this paper differs from the quantity 6 in PS71 except for the 
special case of spheres. This is because PS71 defined t gat to be the damping time for rotation 
about an axis perpendicular to the grain symmetry axis. The subsequent discovery (P79) 
that Barnett relaxation aligns 17 parallel to the symmetry axis of an oblate spheroid implies 
that the definition of t ga , adopted here is appropriate. 

The symmetry of the diffusion coefficients in j x and j y shows that Larmor precession 
causes the angular momentum distribution to be axisymmetric about the the magnetic field 
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as expected (Martin 1971). It might seem that one could exploit this symmetry to reduce 
the number of coupled Langevin equations from three to two in number. For example, one 
might introduce the variable 

jp = sJTl + Jy > 

and then evaluate 

(cos 2 0) = {(jp/j*) 2 ) ( 4 - 9 ) 

by n um erically computing the sample paths of j p and j z . This would amount to a change of 
variables in the Langevin equation. Although are stochastic (hence non-differentiable, 
cf. eq. 2.11) variables, such a change can indeed be performed if one uses “Ito’s formula 
from the stochastic calculus (Gardiner 1990). However, an unfortunate consequence of this 
procedure is that that the resulting Langevin equation for j p contains a term which becomes 
singular when j p = 0. While it seems likely that the Euler-Maruyama method could be 
modified to handle this singularity, it is not clear how much these modifications would 
reduce the potential computational savings — which are in any case no more than 33% 
over the straightforward approach of computing j x , j„, and j z . We therefore adopt the 
straightforward approach here. 

4.2 Benchmark Calculations for Spheres and Thin Disks 

Jones & Spitzer (JS67, see also PS71) obtained an exact solution for the steady-state 
angular momentum distribution of spheres under the influence of sticking gas-grain collisions, 
evaporation, and paramagnetic absorption. The applicability to this paper of the JS67 
solution follows directly from the fact that our diffusion coefficients reduce to those of JS67 for 
spherical grains. Indeed, this correspondence is required by the fact that, for spheres, Barnett 
relaxation does not occur and Larmor-averaging has no effect on the diffusion coefficients for 
collisions and evaporation. In addition, we point out the rather surprising result that the 
JS67 solution is also exact for thin disks. The latter conclusion follows from the fact that R 
depends on the mantle eccentricity only through the ratio Fj./r|| (cf. eq. 4.4-4.7), and this 
ratio is unity for both e m = 0 and e m = 1 (cf. eq. 3.37—3.38). Thus, a thin disk and sphere 
will be aligned with equal efficiency for equal values of the parameters 6 and Td/T g . s 
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5 Of course this does not imply that spheres and thin disks will be aligned equally under a 
given set of physical conditions. The parameter 6 is proportional to the grain volume, V, 
which is zero if e m = 1. Thus the magnetic alignment of an infinitely thin disk is always zero 
and, in this sense, the exact solution for thin disks is only of academic interest. However, 
r x /r|| never differs from unity by more than about 15% (Fig. 2), suggesting that the exact 
solution for spheres and disks may furnish a good approximation for grains of all intermediate 
eccentricities. This possiblity is explored in the numerical calculations presented in Paper 

IL 

Although we will ultimately be interested only in the steady-state angular momentum 
distribution, it is useful to consider the time-dependent distribution, for the purpose 

of estimating various time scales of interest. It is possible to find exact, time-dependent 
solutions which satisfy various initial conditions and which reduce to the JS67 solution as 
f +. oo. Here we will consider an arbitrary but representative example where 

ji{ 0) = 0, i = x,y,z, (4.10) 

i.e., where the initial conditions are deterministic with the grains initially at rest. One can 
easily show that the FP equation is satisfied if the distribution of each component, ji, is a 
statistically independent Gaussian with zero mean value. The time-dependence enters via 
the variances, which are 

„i = = 1 + a - «pi-2(i + ■ < 4U > 

v 2{l + 6) 

for the components normal to B and 

a] _ 1 + (Td/Tg) _ eX p(_2r)} (4.12) 

for the component parallel to B, where r is dimensionless time in units of t gat . Thus the 
distributions of j x and j v relax to a steady state on a time scale of order t gat /{ 1 + 6) while 
the corresponding relaxation time for j z is t gas . The relaxation time is evidently just the 
time scale on which rotational drag causes a grain to “forget” its initial angular momentum, 
this time scale is determined by gas drag alone for the component parallel to B and by gas 
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drag plus paramagnetic absorption for the components perpendicular to B. In the following 
discussion, the term “relaxation time” will refer specifically to t gaa , i.e. to the larger of the 
two relaxation times above. Henceforth we consider only the stationary distribution on the 
assumption that the gas and grain properties are constant over time scales (Table 1). 

In a steady state, the statistic (cos 2 0) on which R depends is given exactly by 

1 * [(^* £ ) 1 ^ sin* 1 " 1 x ^ 2 ~ l] » ifl >° 

(4.13) 

i [(uj )' /I sin-'C-x) 1 ' 1 - l] , ifx<0 

where 

i = a 2 J a] - 1 (4.14) 

(cf. PS71, eq. 23). The steady state autocorrelation function of ;',(r) can be shown to be 

(Ji( T )ji( T> )) = a i eX P H r - T '\/ T cor,i] , ( 4 i5 ) 


where the dimensionless correlation times in units of t ga , are given by 


'Tcor,! — \ 


t h if * = x >» 


(4.16) 


1 


if i = z. 


Henceforth, we will take the term “correlation time,” to mean the smaller of these time 
scales, i.e. 

Tear = ( 4 - 17 ) 


l+<' 


We expect, rougly speaking, that the minimum time scale for fluctuations in cos 2 /? will be 
of order r^. 


Numerical estimates of (cos 2 0) for spheres and disks were obtained for comparison 
with exact solutions as follows. Given values of the parameters 6 and Td/T g , we calculated 
N = 50 values of (cos 2 0) from N identical trial calculations. In each trial, we took the initial 
conditions as in eq. (4.10) and then used the Euler- Maruyama scheme, eq. (2.14), to generate 
a sample path of j (t). The Box-Muller algorithm was used to generate the Wiener increments 
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from a sequence of uniform pseudo-random numbers and the latter were obtained, in turn, 
using the portable random number generator referred to as “RAN3” in Press et al. (1986). 
In each trial, we integrated the sample path initially for 100 relaxation times to let the grain 
forget its initial conditions. Thereafter the time averaging of cos 2 0 commenced and the 
Langevin equation was integrated over an additional, dimensionless “averaging time,” r avg . 
The error in each trial integration was characterized by E re i , the relative error in (cos 2 0) 
determined by compairing the trial result with the exact solution. We expect the accuracy 
of (cos 2 0) to depend on r avg , on the dimensionless time step, A r, and on systematic errors 
caused by various effects such as spurious correlations in the random number generator. We 
now consider each of these effects individually. 

The dependence of errors on the averaging time is illustrated in Figure 5, where we have 
plotted the root-mean-square value of the N = 50 values of E re i as a function of r avg . All of 
the calculations in Fig. 5 pertain to spheres or thin disks with Td/T g = 0.5, 6 — 1, and used 
a time step At = 1 x lO' 3 ^. The straight line in Fig. 5 is a least squares fit to the plotted 
points. Its slope, -0.47, is consistent to within the statistical errors in E rm , with the scaling 
law Em* oc (Tavj/Tcor) -1 ^ 2 . Of course this is just the scaling one expects on the basis of the 
Central Limit Theorem: the value of (cos 2 0) at the end of a trial is effectively the mean of 
M statistically independent estimates, where each independent estimate corresponds to the 
average over one correlation time and M = T avg /r cor is the number of correlation times in 
one trial integration. 

The dependence of E rm . on the time step is shown in Fig. 6. Each rms error was 
determined as in Fig. 5 from the results of 50 identical trial integrations. The total integration 
time for all points in Fig. 6 was r avg = 5 x 10 4 r cor and the physical parameters were chosen to 
be 6 = 1, Td/T g = 0.5. For large At, where the accuracy in (cos 2 0) is limited primarily by 
the time step, E rmt scales approximately as ( At/t co , ) -m , where m « 1. As At decreases, the 
rms error approaches a constant, reflecting the fact that the accuracy in (cos 2 0) is eventually 
limited by the integration time. Thus the minimum useful time step is about At = 10' 2 for 
an integration time of r avg = 5 x 10 4 r COT . Note that the slight upturn in E rm . for the smallest 
time step, At = 1 x I0 -3 T c<>r , is not statistically significant: The increase in log Ermt-, about 
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0.05, is approximately equal to l/>/3V, where N = 50 is the number of independent trials 
used to deter min e E rmt . That is, the upturn is only a “1<t fluctuation. 

Systematic errors due to spurious correlations in the random numbers, etc, should reveal 
themselves as residual mean errors in calculations where the two effects discussed above are 
small. An upper limit on any systematic errors which may be present in our numerical results 
can be inferred from the data in Figures 7 and 8. The histograms depict the distributions in 
Erei for N = 100 identical trials with Ar/r^ = l.Ox 10 -3 and two different integration times. 
All of the calculations in Fig. 7 and 8 pertain to the case S = 3 and Td/T g = 0.9. In Fig. 7, 
where r avg = 6 x lO 4 ^, the relative errors have mean value ~E rt i = -1 x 10 -3 and standard 
deviation a E = 5 x 10~ 3 . In Fig. 8, where r avg = 1 x lO 6 ^, the corresponding statistics 
are ~E re i = -6 x 10“ 6 and a E = 1 x 10 -3 . We conclude from the smallness of |f? re /| l°E that 
any systematic errors which may be present are much smaller than the statistical errors. 
The “production calculations” presented in Paper II use a time step of At = 1 x 10 3 r cor 
and an integration time of r avg = 6 x 10 4 r eor 6 . The histogram in Fig. 7 therefore reflects 

6 The reader may notice that our integration time does not warrant such a small time step 
(cf. Fig. 6). We choose a small time step to allow for the possibility of averaging our results 
with additional integrations if computer time becomes more plentiful in the future. 
the distribution of relative errors in our production runs. We conclude that our production 
calculations have a 2 a error in {cos 2 ft) of about 2%, which translates to a 2cr error in R of 
about 6% when R = 0.1. When R <C 0.1, the relative error in R is generally much larger 
than 6% because (cos 2 /?) * 1/3 (cf. eq. [1.1]). This effect is shown in Fig. 9, where the 
symbols depict il- values for spherical grains determined by numerical integration and the 
solid curves give the exact analytic solutions for comparison. The rms error in R (not log R) 
for all points in Fig. 9 is 1.3% but the point with R = 5 x 10~ 3 has an error of about 25%. 
We could reduce the errors in very small R values somewhat by integrating longer and with 
a smaller time step, but there does not appear to be any incentive for doing so at present. 
The 25% error in R = 5 x 10 -3 is an extreme example which, in any case, has no observable 
consequences. For example, the linear polarization in fax-infrared or submillimeter emission 
from aligned, oblate spheroids made of astronomical silicates (Draine 1987) would be about 
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60#% for an axis ratio b m /a m = 2 if we neglected (unrealistically) effects which reduce the 
effective degree of alignment, such magnetic field line inhomogeneities along the line of sight. 
Even in this optimistic picture, the polarization corresponding to R = 5 x 10“ 3 would be 
near the limit of detectability (e.g., Hildebrand et al. 1992) and a 25% error in R would be 
unobservable in this case. 

5 SUMMARY 

1. We have demonstrated a new method for obtaining accurate predictions of grain align- 
ment theories which uses numerical integration of the Langevin equation. The new 
method is easy to implement and superior in several respects to Monte Carlo tech- 
niques or the direct solution of the Fokker-Planck equation. 

2. We have calculated the Fokker-Planck diffusion coefficients for physical processes rele- 
vant to super-paramagnetic grain alignment in molecular clouds, on the ad hoc hypoth- 
esis that the grains do not rotate suprathermally. These calculations will be applied 
in the next paper of this series to derive observational tests of this hypothesis. 

3. We have developed a quantitative theory for the effects of thermal emission and far- 
infrared absorption on the alignment of nonspherical, dielectric grains. These processes 
generally have a negligible effect on grain alignment in molecular clouds. 

4. We have performed benchmark calculations to calibrate the errors in our computational 
methods. We have determined the dependence of errors in the Rayleigh reduction factor 
on parameters in our numerical integration scheme. Any systematic errors that may 
be present in our calculations correspond to spurious predictions of nonzero alignment 
which are below the present limits of detectability. 
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APPENDIX A 

STICKING COLLISIONS ON AN INFINITESIMAL SURFACE 

PATCH 

We compute the incremental contributions to the mean torque and diffusion tensor 
which are caused by sticking collisions onto an infinitesimal surface patch of the grain. The 
results of this appendix apply to grains with arbitrary shapes, angular velocities, and gas- 
grain drift velocities. (The effects of gas-grain drift are included for reference in future 
work.) The incremental contributions are integrated over the surface of an oblate spheroid 
in Appendix C to obtain the total diffusion coefficients for sticking collisions. 


A.l Differential Flux of Colliding Molecules 

Let v = Vth 8 be the velocity of a molecule in the gas frame (i.e., its thermal velocity), 

iatiImq an t \ 

1/2 


where 8 is dimensionless and 


Vth 


■(?) 


(Al) 


is the gas thermal speed. If the distribution of thermal velocities is Maxwellian at the gas 
temperature, then the mean flux of colliding particles with velocities in the range («, 8 + <Ps) 
is nvth f c (a)<p8, where 


/c(s) = { 


7 r -3 / 2 (s p — s)*n exp (— a 2 ) if (s — s p )-n < 0 

0 otherwise 


(A2) 


and 

8 P = (tij -I- f1Xr)/v t K M 3 ) 

is the total tr ans lational plus rotational velocity of the patch in the gas frame in units of 


A. 2 Incremental Mean Torque 

Consider an infinitesimal surface patch of the grain with area dA, unit outward normal 
n, and position vector r relative to the grain CM. When a molecule with velocity v t h 8 sticks 


36 


to the patch, the grain internal angular momentum increases by 


SJ{a) = mvth [rx(s - «,*)] , (^ 4 ) 

where v d = v lh a d is the CM velocity of the grain in the gas frame. Multiplying expression 
(A4) by the mean flux of colliding particles and then integrating over all thermal velocities 
gives the patch's contribution to the mean torque due to sticking collisions. Its cartesian 
components are 

d (A Ji) c = n v th dA j 6Ji{a) f c (a) <Pa, ( A5) 

where in this and subsequent expressions integrals over a are assumed to be definite integrals 
over all of s-space. The integral in eq. (A5) can be evaluated analytically to give 

d(AJ,) c = nmv ( \ [G(a p *n) (r Xn),- - F(a p 'h) (r Xa d )i ] dA, (A6) 

where the functions F and G are defined below in §A.4. 

A.3 Incremental Diffusion Tensor 

The patch’s incremental contribution to the diffusion tensor has cartesian components 


d (A Ji A Jj) c = n v th dA J 6Ji(a) 6Jj(a) f c (a) <Pa. (A7) 

The integral can be evaluated analytically by straighforward calculation to show that 

i (AJ. A J,), = £<i( AJ< A ■/,)<*> , ( A8) 

k=l 

where 

d (A Ji A J,)™ = i nm 2 v^F(a^n)T^dA, ( A9) 

d{AJi AJj)™ = nm 2 Vf h {H(a p -h) - ^(sp-n)] T^dA, (A10) 

d (A Ji AJj)f ] = -nmV/A V^f^A, (All) 
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d (A Ji A = nm J 4F(^,-n)5-J 4) d/l. ( A12) 

The function H is defined below and T (1) , . . . , T (4) are second-rank tensors with components 

(^13) 

= (rXh)i(rXh)j, (- 414 ) 

= (r Xn),(rX«j)j + (rXJjMrXn),-, (A15) 

and 

T^ ] = (r X «d),-(r X 8 d )j, (^16) 

where Sij is the Kronecker delta symbol. 


A.4 Special Functions 

The function 


f <*> s 2 


x(l +erfx) + -^=exp(-x 2 )j (^17) 

is just the dimensionless particle flux integrated over all thermal velocities for a patch whose 
velocity in the gas frame has a component x = « p *n along its own outward normal. The 
limiting behavior of F for large find small |x| is 

^x- 2 exp(-i 2 ) a 


— oo 


F(x) ► < 27 ? + 2 X 


{ X + 47? X 2e xp(-x 2 ) X-^+OO. 


The function 


G(z) 3 -- 


1 + 7 P &» ! ) 




is defined in terms of the incomplete gamma function, 

P(a,u) = f t a _1 exp(— t) dt, 

I ^aj 70 


(A18) 


(A19) 


(420' 
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and has the limiting forms 


G(x ) — 


The function 

H(x) = ^ x 

has the limiting behavior 


H(x) 




X — ► — oo 


.1 £_ 

4 


x — * 0 


“2 + S77 ex P(— x2 ) x-^+oo. 



+ 27* (1 + x 2 )exp(-x 2 ) 


47? ex P(~ 


* 2 ) 


x — ► —00 


< ife + J* *-° 

\ x + 577 exp (-x 2 ) x-»+oo. 


{A2\) 


(A22) 


(A23) 
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APPENDIX B 

THERMAL EVAPORATION FROM AN INFINITESIMAL SURFACE 

PATCH 

Here we compute the contribution of an infinitesimal surface patch to the diffusion coef- 
ficients for thermal evaporation. The results of this appendix apply to grains with arbitrary 
shapes, angular velocities, and gas-grain drift velocities. The total diffusion coefficients for 
thermal evaporation from an oblate spheroid are derived in Appendix C by integrating the 
results of this appendix over the grain surface. 


B.l Differential Flux of Evaporating Molecules 


Let v' = v ev a' be the velocity of an evaporating molecule in the frame of reference of 
the surface patch, where s' is dimensionless and 



(B 1 ) 


is a characteristic thermal speed for evaporation. Consider first the hypothetical situation 
where the gas and grain temperatures are equal and the patch is at rest in the gas frame. 
In this special case, the principle of detailed balancing requires that every evaporation with 
velocity v' should be accompanied, on average, by a sticking collision with velocity —v . The 
velocity distribution of evaporating molecules is therefore proportional to 


fev (S') = { 


| a' -ft exp [- (s') 2 ] 

0 


if «'*n > 0 
otherwise, 


m 


which is equal, apart from a normalization factor discussed below, to the flux of sticking 
collisions, eq. (A2), with a' = -a and the patch velocity, * p , set to zero. Although it was 
derived assuming that T d = T g , expression (B2) is correct for all T d and T g because the evap- 
oration spectrum is obviously independent of the gas temperature. We will also assume that 
expression (B2) is correct for a rotating surface patch. The neglect of any possible dynamical 
coupling between grain rotation and evaporation should be an excellent approximation for 
grains rotating at thermal angular velocities, because the centrifugal potential of a molecule 
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rotating along with the grain surface is only of order ( m/M* ) kT g . We have defined /«, to 
be normalized so that 

y.u»')<fv=i, (S3) 

where in this and subsequent expressions it is assumed that integration over a 1 means in- 
tegration over all of s^space. The differential flux is therefore $ f ev (s 7 ), where $ is the 
mean n um ber of evaporations cm -2 s -1 integrated over all a' . If the evaporation rate is 
independent of position on the grain surface, which is the case of interest in this paper where 
suprathermal rotation is not considered, then $ is set by the requirement that the total 
collision and evaporation rates integrated over the surface must balance. Thus 

| _ nv„ ffWn) 4A (B4) 

iS* 

where v th is the gas thermal speed, S is the total grain surface area, the integral is over the 
grain surface, and the function F(sp'h) is defined in Appendix A. 

B.2 Incremental Mean Torque 

A molecule which evaporates with dimensionless velocity s' relative to the rotating, 
drifting grain surface has a (dimensional) velocity 

t) = u ro «' + f?Xr (£5) 

relative to the grain CM and its departure changes the grain internal angular momentum by 

8J{J) = -m { u ev (rxi')+rx (f?Xr) } . (B 6) 

The incremental m ean torque due to evaporation from a single surface patch therefore has 
cartesian components 

d < A Ji) m = *dA J 6Ji(a') /,„(«') da 1 . (£7) 

Straightforward evaluation of the integral above shows that 

d ( A = — m $ | ~^jr Vev ( r * ”)» + dA > ( B8 ) 

where denotes the scalar product of fi with the second-rank tensor (cf. Ap- 

pendix A). 
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B.3 Incremental Diffusion Tensor 

The patch’s incremental contribution to the diffusion tensor has cartesian components 

d (A Ji A J,) tv = $dA J SJi(s') 8J,(a') f ev (a') da'. (BO) 

When the integration over a ' is performed, the cartesian components are found to be 

d ( A J t A Jj)" = E d (A Ji A Jj)™ , (BIO) 


k= 1 


where 


[rV+r^iA, 


(B 11 ) 


d(AJ, AJ,)™ = ^nSv.A 


(rXn)j (TW-n). + (T( 1 >-I7).(rxn) j 


and 


d(&Ji &J,)M = m 1 i (T(')-n) . 


<M, (512) 


(513) 
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APPENDIX C 

DIFFUSION COEFFICIENTS FOR STICKING COLLISIONS AND 
EVAPORATION FROM AN OBLATE SPHEROID 


C.l Oblate Spheroidal Coordinates 

Here we compute the diffusion coefficients for sticking collisions and thermal evaporation 
by integrating the incremental diffusion coefficients from Appendixes A-B over the surface 
of an oblate spheroid. In the present paper we will assume that the grain center of mass 
is at rest in the gas frame ( v d = 0) and that the grain is rotating about its symmetry axis 
due to the effects of internal dissipation. The surface integrals are easy to evaluate if one 
specifies the position of a surface element by its oblate spheroidal coordinates, (*/,£)> where 
rj € [-tc/ 2 ,it/ 2] and $ € (0,2*] (see, e.g., Arfken 1970). Thus the vector r(»?,0 from the 
grain CM to the mantle surface at (77, £) has body-frame coordinates 


The length of r is 


x b (rj,{) = b m cos*/ cos £ 
y b (v,0 = cos rj sin £ 
z b (r},Z) — ° m sini/. 

1 1/2 

r{v, i) = [cos 2 rj + (1-4) sin2 »?| . 


the unit outward normal to the surface at r hats body-frame coordinates 

n£(? 7, 0 = [sin 2 77 + (1 - 4) cos 2 17] (1 - 4) 1/2 cos 7 cos 4 


n*(»/,0 = [sin 2 T) + (l — 4) cos 2 77 
n‘(i7,{) = [sin 2 77 + (1 — 4) coa2 V 

and the differential element of surface area is 

[(i- 4 ) 


-1/2 

-1/2 


(1 - 4) 1/2 coST ? sin £ 


smt7, 


— tfn 4 


+ sin 77 


m 


1/2 


cos r 7 dr} d£. 


(Cl) 


(C2) 


(C 3) 


(C 4) 
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C.2 Sticking Collisions 

C.2. 1 Mean Torque 

The mean torque is obtained by using expressions (Cl)-(C4) to write the incremental 
diffusion tensor of a surface patch, eq. (A6), as a function of its oblate spheroidal coordinates 
and then integrating over 17 and £. First we note that 8 p -h = 0 everywhere on the mantle 
surface because v d = 0 and the grain is a surface of revolution which is rotating about its 
symmetry axis. Noting also that the second term in (A6) vanishes when v d = 0, and that 
G(8 p 'h) = G(0) = -1/4 is just a constant, we see that the mean torque due to collisions has 

cartesian, body-frame components 

(&Ji) c = ~^nmv 2 tll J (r Xn)J dA. (£5) 

The integral above vanishes for any surface of revolution. We conclude in agreement with 
PS71 that the mean torque due to sticking collisions is zero, 

(a j?) c = 0 , m 

for any surface of revolution which is rotating about its symmetry axis. 

C.2. 2 Diffusion Tensor 

To integrate expression (A8) over the grain surface, we first make use of the fact that 


s d = 0 and s p *n = 0 as noted above, from which it follows that two of the tensors in (A8) 

vanish, 


J(3) = 7 ( 4 ) = 0, 

{Cl) 

and also that the functions 


F(vn) = 575 

(C 8) 

and 


"<V*> = 

(C9) 


are just constants. The diffusion tensor for sticking collisions therefore has body-frame 


components given by the surface integral 

(AJ*AJ‘) c = ^i=nm 1 v,\ J (CIO) 
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The integrals can be evaluated analytically for each i and j by first writing out the tensor 
components as functions of (tj, if) using eq. (Cl) and (C3). Straightforward calculation shows 
that (A^AJ/) e is diagonal with components 

((AJ{) 3 ) c = ^nm 2 6^vf fc r||(e m ) (Cll) 

and 




(C 12 ) 


where the subscripts “||” and “X,” respectively, denote body frame axes parallel and normal 
to the symmetry axis. The shape dependence of the diffusion tensor is expressed in the 
functions 


r,|M = i { 3 + 4(1 - 4)sM - e- 1 [l - (1 - O’jfe,,)]} , (C13) 


and 


ri(em) = 4 { 7 - + (> - «!)»(«-) + 0 - 2 4) [l + [l - (1 - e™) J 9(«™)]]} • 

(C14) 

where 

(£15) 




C.3 Thermal Evaporation 

C.3.1 Mean Torque 

The integration of the incremental mean torque over the mantle surface is simplified by 
noting that, because s p *n = 0 everywhere on the surface, the mean evaporation flux is just 

« = (£16) 

20T 

The first term in eq. (B8) vanishes when integrated over the mantle surface as noted above. 
Integration of the remaining term shows that 

= 0 

Wi)„ = o < cu ) 

(A Jx)^ = ||(e m )£I, 
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where we used the fact that 17 = fli* due to the effects of internal dissipation. Thus 

(A J)„ = -J/V-, < C18 > 

where the gas damping time is defined to be 

t = hi (C19) 

9at 4 y/fr nmu t /,6^r||(e m ) 

and we have dropped the superscript “b” in eq. (C18) because it is written in a form which 
is valid in any coordinate system. 

C.3.2 Diffusion Tensor 

Our assumption that the grain rotational energy is of order kT g implies that magnitudes 
of expressions (Bll), (B12), and (B13) are in the approximate ratio hyJm/M d :m/M d . Thus 
the diffusion tensor for evaporation is approximately 

= -^=nm 2 v th vl v J + T/f' b) {T},^ dA{ri,Z), (C20) 

where, in neglecting the contributions due to expressions (B12) and (B13), we have intro- 
duced a relative error which is typically of order 10 -4 for large grains with 6 m ~0.1 fim. 
Comparing eq. (C20) with eq. (CIO), we see that the diffusion tensors for collisions and 

evaporation axe related by 

(AJ‘A7‘) tv =|(A7‘A7‘) t . «* U 

This result also could have been deduced from the principle of detailed balance (see, e.g., 
JS67), i.e. without performing the integration in eq. (C20). We derived expression (C21) 
directly as a check on the complicated series of calculations which yielded the diffusion tensor 
for collisions. 
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Table 1 


Time Scales for an Oblate Spheroidal Grain® 


Symbol 

Definition 

Value 

trot 

Grain Rotation Period 

6.3xlO- s T-y 2 o^,l 5 iiS,,.,eSJ 1 

tlar 

Larmor Precession Period 

3.7 x 10 s Bll (x'-s)’ 1 *>m,- 5 ©« 

tbar 

Barnett Relaxation Time 

3.3 x 10 7 KZ \ 3 T~\ am, -s ©L ^ (* ~ l)' 1 

tgas 

Gas Damping Time 

6.2 x 10 9 n 4 " 1 r;} /, fl Wl - 5 ©« 


“Time scales in seconds for an oblate, spheroidal, core-mantle grain rotating about its sym- 
metry axis with kinetic energy kT g / 2. The notation “5-s” means B/ 10 5 G and similarly 
for other quantities in Gaussian cgs units. Core and mantle densities are assumed to be 
p c = 2.5 g cm -3 and p m = /> c /3. Magnetic properties x' and K are volume averages over 
core plus mantle and are scaled appropriately for ordinary (non-superparamagnetic) sub- 
stances. For super-paramagnetic grains, ti ar and U aT should be reduced by a factor £ 10 . 
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FIGURE CAPTIONS 


Fig. 1 — Relationship betwee the basis vectors (x, y, z) of the gas frame and (x\ y b , z b ) of 
the body frame. The transformation of vectors and tensors between the two frames depends 
on the Eulerian angles, /?, <f>, and rj>. 

Fig. 2 — Eccentricity factors, T ± and T||, defined by eq. (3.34)-(3.35), plotted vs. the mantle 
eccentricity. The efficiency of grain alignment by collisions, evaporation, plus paramagnetic 
absorption (§4) depends on grain shape only via the ratio T x /U || , which is also plotted. 

Fig. 3 — Photoabsorption damping time, eq. (3.58), for silicate grains immersed in a thermal 
radiation bath with temperature T ra d- For simplicity we consider bare silicate grains with 
density p e = 2.5 g cm -3 and axis ratios of 1:1 (solid), 2:1 (long dash), and 5:1 (short dash). 
In each case we consider three radii, b c = 0.02, 0.1, and 0.5 pm. The largest and smallest 
damping time for each axis ratio correspond, respectively to the largest and smallest radius. 

Fig. 4 — Similar to Fig. 3, but the quantity plotted is the photoemission damping time, 
eq. (3.67). 

Fig. 5 — Symbols: Root- mean- square value of the relative errors in (cos 2 (3) obtained from 
50 identical trial calculations for spherical grains. E rm * is plotted as a function of T avg / t C0T , 
the averaging time in units of the correlation time. All calculations pertain to spheres with 
T d /Tg = 0.5, 6 = 1 and used a dimensionless time step of At/t cot = 1 x 10 -3 . Solid curve: 
Least squares fit of a straight line to the symbols. The line has slope —0.47. 

Fig. 6 — Similar to Fig. 5, but E rmt is plotted as a function of the integration step size in 
units of the correlation time. All other parameters have values as in Fig. 5, except that the 
averaging time is 5 x lOV^. 

Fig. 7 — Histogram: Distribution of the relative errors in (cos 2 /?) for 100 trials with 
A t/tcot = 1 x 10~ 3 and r OVJ /r COT = 6 x 10 4 . All calculations pertain to spheres with 



Td/Tg = 0.9 and 6 = 3. The plotted error distribution has mean value E rm , — 1.1 x 10 

and standard deviation = 5.1 x 10' 3 . Dashed curve: Gaussian distribution with zero 
mean, standard deviation equal to <r E , and normalization chosen to have the same area as 

the histogram. 

Fig. 8 — Similar to Fig. 7, but the averaging time for each trial is T avg /T C0T = 1 x 10 6 . 
The error distribution indicated by the histogram has mean value ~E rm . = -5.8 x 10 6 and 
standard deviation &e = 1.3 x 10 3 . 

Fig. 9 — Rayleigh reduction factor for spherical grains plotted vs. S. Symbols: Numerical 
results obtained by integrating the Langevin equation with At / r cor = 1 x 10 -3 and T axjg /T cor = 
5 x 10 4 . Solid curves: Exact solution for spheres obtained from eq. (4.14). The root-mean 
square relative error in the numerical values of R (not log R) is 1.3%. Notice that the relative 
errors are larger for small R- values because (cos 2 /?) -► 1/3 as R -* 0 (cf. eq. [!.!])• 
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